The Elastic Behavior of Entropic "Fisherman's Net" 
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A new formalism is used for a Monte Carlo determination of the elastic constants of a two- 
dimensional net of fixed connectivity. The net is composed of point-like atoms each of which is 
tethered to six neighbors by a bond limiting the distance between them to a certain maximal 
separation, but having zero energy at all smaller lengths. We measure the elastic constants for 
many values of the ratio 7 between the maximal and actual extensions of the net. When the net is 
very stretched (7 ~ 1), a simple transformation maps the system into a triangular hard disks solid, 
and we show that the elastic properties of both systems, coincide. We also show that the crossover 
to a Gaussian elastic behavior, expected for the non-stressed net, occurs when the net is more loose 
(7 ~ 3). 
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Materials like rubber and gels are formed when poly- 
mers or monomers are cross-linked into macroscopically 
large networks. Due to the small energetic differences (of 
the order of kT) between the allowed microscopic con- 
figurations of these materials, their physics is primarily 
determined by entropy, rather than energy. This has been 
recognized long ago, and the peculiar physical properties 
of rubber and gels, in particular their great flexibility, 
are attributed to this microscopic feature. The classical 
theories of rubber elasticity, for instance, deal with Gaus- 
sian networks in which the internal elastic energy is com- 
pletely ignored and the strands between cross-links are 
viewed as entropic springs These theories, however, 
do not explain well the elastic behavior of networks of 
certain types. Perhaps the most known unresolved prob- 
lem in this field of research, is the question of the criti- 
cal elastic behavior of random systems near connectivity 
threshold. Most of the numerical works which aimed to 
investigate this issue during the last twenty years, con- 
cerned with the energetic elasticity j^. Recent studies 
however, suggested that close to the gel-point elas- 
ticity is dominated by its entropic component. A com- 
pletely different aspect of entropic elasticity which has 
been studied much less, is the behavior of highly con- 
nected networks, well above their connectivity threshold. 
The classical theories are inappropriate in this connectiv- 
ity regime, since the strands between the junctions are 
very short and do not resemble Gaussian springs. 

When the boundary of a thermodynamic system is ho- 
mogeneously deformed, the distance between any two 
boundary points which prior to the deformation were sep- 
arated by R, becomes 



1/2 



(1) 



where the subscripts denote Cartesian coordinates and 
summation over repeated indices is implied. The quan- 
tities r]ij are the components of the Lagrangian strain 
tensor, while Sij is the Kronecker delta. The elastic be- 
havior of the system is characterized by the stress tensor, 
aij, and the tensor of elastic constants, C'ijki, which are 



the coefficients of the free energy density expansion in 
the strain variables 

/({'7}) = /({O}) + cTijVij + ^CijkiVtjVki + . . . . (2) 

Measuring the elastic constants is much more difficult 
in entropy-dominated systems than in energy-dominated 
ones. In the latter one needs to calculate energy varia- 
tions around well defined ground states. In the former, 
on the other hand, different microscopic configurations 
possess similar energies. Entropy in this case is essen- 
tially the (logarithm of the) number of allowed micro- 
scopic configurations. Measuring the variations of this 
quantity in response to external deformations applied on 
the system, is usually very complicated. In order to sim- 
plify this task, and due to the fact that the exact energy 
details are quite irrelevant in entropy-dominated systems, 
the inter-atomic interactions in such systems are often 
modeled by "hard" potentials. Excluded volume effects, 
for instance, can be modeled by the hard spheres repul- 
sion, while chemical bonds can be replaced by inextensi- 
ble ( "tether" ) bonds which limit the distance between the 
bonded monomers, but have zero energy at all permitted 
distances The energy of all the microscopic config- 
urations in such models, which are called "athermal" , is 
the same, and their physics, therefore, is exclusively gov- 
erned by entropy considerations. It is interesting to note 
that although athermal models have been investigated 
quite extensively in polymer and soft matter physics, the 
elastic properties of many of them are not well under- 
stood. Hard spheres systems, for instance, are studied 
for already more than 40 years They were, in fact, 
the first systems for which Metropolis et al. performed 
the first Monte Carlo (MC) simulations in 1953 §. The 
phase diagram of hard spheres, which is a function of a 
single parameter, their volume fraction, had been fully 
explored both in simulations and experiments Q. Yet, 
despite of the numerous works dedicated to elasticity of 
hard spheres ||^ , the accuracy of the values of their elastic 
constants still leaves much to be desired. 

In the canonical ensemble, the elastic constants can be 
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related to the mean squared thermal fluctuations of the 
stress tensor components (just as the heat capacity is pro- 
portional to the mean squared energy fluctuation). This 
relation, first expressed by Squire et al. can be used 
for a Monte Carlo determination of the elastic constants. 
The method is known as the "fluctuation method" . The 
instantaneous stress, measured at a given microscopic 
configuration, is associated with the mean force (aver- 
aged over the entire volume) acting on the atoms [ p^ . 
The local forces originate from external potentials and 
inter-particle interactions. In entropy-dominated sys- 
tems, these forces are usually very small. They become 
extremely large only over very short time intervals when 
atoms come to the close vicinity of each other or when 
bonds are sufficiently stretched. Model with hard poten- 
tials can be regarded as the limiting case in which these 
time intervals vanish, while at the same time the instan- 
taneous forces become infinitely large, keeping the rate of 
momentum exchange between atoms fixed. It is obvious 
that the stress in such systems must be related to the two- 
point probability densities of contact between spheres 
and occurrence of bond stretching, while the elastic con- 
stants (stress fiuctuations) must be related to the cor- 
responding four-point probability densities. Indeed, we 
have recently succeeded to formulate the exact relations. 
We obtained expressions enabling a direct measurement 
of the entropic contribution to the elastic constants, and 
demonstrated the accuracy and efficiency of the method 
using this formalism on three-dimensional hard spheres 
systems [pT| . In this paper we apply this new formalism 
to measure, by means of MC simulations, the stress and 
elastic constants of topologically simple networks. We 
consider a "toy model" consisting of a two-dimensional 
(2D) network of atoms forming a triangular "fisherman's 
net" (FN): atoms are point-like, i.e., have no excluded 
volume, and each one of them is connected to six neigh- 
bors by a "tether" limiting the maximal distance between 
the atoms, but otherwise not exerting any force. The FN 
is a highly connected network, whose physical behavior 
is entropy-dominated. Very few studies were devoted to 
systems of this type, and it is indeed quite clear that the 
determination of the elastic constants of systems similar 
to the FN, is far from being trivial. 

The FN is six-fold symmetric when it is equally 
stretched along all the spatial directions. Its elastic prop- 
erties in this reference state should be as of an isotropic 
system Its stress tensor is diagonal with the ele- 

ments a — a XX — — —P, where P is the negative 
external pressure (stretching) one needs to apply to the 
boundaries in order to balance the forces exerted by the 
net. Only four elastic constants of the net do not van- 
ish. Cll — Cxxxx^ C22 — ^yyyyi — Cxxyy — Cyyxx 
and C44 — Cxyxy — C'yxyx — ^xyyx — ^yxxy- DuC tO the 

isotropic nature of the system, only two of them are inde- 
pendent, and they satisfy the relations: Cn = C22, and, 
2C44 = Cll — C12 It is quite common to describe 

the elastic properties of isotropic systems in terms of the 
shear modulus, /x, and the hulk modulus, k, defined by: 



^ = C44 — P, and, K = 5(6*11 -I- C12). When these quan- 
tities are positive, the isotropic system is mechanically 
stable 0. 
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FIG. 1. Configurations corresponding to different values of 
the ratio 7 between the maximal and actual extensions of the 
net: (a) 7 = 1.1, (b) 7 = 1.5. Only part of the net is shown 
in the figures. 

Our simulations were performed on systems consisting 
of 1600 atoms which were bonded to form a triangular 2D 
net. The topology of the net is such that the mean posi- 
tions of the atoms form a regular triangular lattice with 
lattice spacing 60 j while each pair of nearest neighbor 
atoms is connected by a tether whose maximal extension 
is 6 > 60. Periodic boundary conditions, which fixed the 
volume and prevented the net from collapsing, were ap- 
plied. We denote by 7 = the ratio between between 
the maximal and actual extensions of the net. Typical 
equilibrium configurations corresponding to two values of 
7 are depicted in Fig. |^. We generated the MC configu- 
rations using a new updating scheme, recently proposed 
by Jaster , in which the conventional Metropolis step 
of a single particle is replaced by a collective step of chain 
of particles. At each MC time unit we made 1600 move 
attempts (with acceptance probability ^ 0.7), where at 
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each attempt a new atom was selected randomly. (On the 
average, each atom was chosen once in a MC time unit.) 
Correlations between subsequent configurations were es- 
timated from the autocorrelation function of the ampli- 
tude of the longest-wavelength phonon in the systems 
(both longitudinal and transverse phonons were checked) . 
For all 7 values, we found that after less than 1000 MC 
time units, the memory of the initial configuration is com- 
pletely lost. We measured the stress and elastic constants 
for many values of 7. For each 7, we averaged the relevant 
quantities over a set of 1.5 x lO'' configurations separated 
from each other by 3 MC time units. We also evaluated 
the standard deviations of the averages. The error bars 
appearing in the graphs which present our results corre- 
spond to one standard deviation. More technical details 
of the simulations, as well as a detailed explanation on 
the formalism used in this work, were given in another 
publication 

When the net is fully extended (7 = 1), atoms cannot 
leave their mean lattice positions. Entropy, therefore, 
vanishes, while the stress and elastic constants diverge. 
For slightly larger values of 7, atoms are restricted to 
small thermal fiuctuations around their lattice positions, 
as in Fig. [l] (a) . A similar atomic-level picture appears in 
hard disks (2D "hard spheres") solids for densities prox- 
imal to the close-packing density. In fact, the FN and 
the hard disks (HD) problems are closely related: In the 
latter (HD) the centers of the disks are not allowed to 
approach their neighbors a distance smaller than a, the 
diameter of the disks, while in the former (FN) atoms 
are not allowed to depart from their neighbors a distance 
larger the maximal extension of the bond, b. For HD 
solids, one can define the ratio S = a/bo < 1 between the 
diameter of the disks, a, and the mean lattice separation, 
60 ■ In the limits 7^1 and 5 — > 1 (corresponding to the 
FN and HD problems, respectively), the elastic constants 
of both systems coincide, as can be seen from the fol- 
lowing argument: Let HpN and Hhd be phase spaces of 
allowed configurations of a FN with a certain value of 7 
and of a HD solid with 6 = I/7, respectively. Each con- 
figuration in one of these phase spaces can be described 
by the set {u^} of deviations of either the atoms of the 
net or the centers of the disks from their mean lattice 
positions. In the 6 ^ 1 asymptotic regimes, we can 
assume that the size of all the deviations is much smaller 
than the lattice spacing, 69. One can easily check that 
if the set {u^} represents an allowed microscopic con- 
figuration of the FN, then the set {— u^} almost always 
corresponds to an allowed configuration of the HD sys- 
tem. Moreover, by this transformation we can generate 
almost all the configurations of Hjjd- The measure of 
the subgroup of configurations for which the mapping 
{u;} < — > {— Ui} between the two problems does not ap- 
ply, diminishes proportionally to {ui'^)/bo^. Thus, the 
mapping {u.^} < — > {— u^} is asymptotically a one-to-one 
transformation from HpN onto Hhd- Since for both sys- 
tems the Helmholtz free energy F is equal to — fcTln |H|, 
where |H| is the volume of the 2iV-dimensional config- 



uration phase space {N is the number of atoms), and 
since the Jacobian of the above transformation is unity, 
we readily find that the free energies i^HD and ^fn of the 
HD and FN systems, respectively, are related by 



FFN(Af,7) -^HD(Af,'5 = 1/7), for 7^1. 



Suppose now that both systems are slightly deformed 
from their reference states. The displacements of the 
atoms from their mean lattice positions can be divided 
into the set {u^} of thermal fluctuations and the set 
{vi} of small changes in mean lattice positions caused 
by the deformation. The transformation between Hfn 
and Hhd, in this case, maps both {uj} to {— Uj} and 
{vi} to {— Vi}. The {vi} mapping is equivalent to the 
reversal of the strain applied on the system. We therefore 
find that i^HD and FpN will be equally modified, provided 
that opposite strains are applied on the FN and HD sys- 
tems. The following asymptotic relations between the 
stress and elastic constants of these systems follow im- 
mediately: crFN(7) - -Phd(1/7), kfn(7) - khd(1/7) and 
C'fn44(7) — C'hd44(1/7)- These relations are very use- 
ful since the asymptotic expressions for Phd, '^hd and 
C'hd44 are available p6| , and can be used to find the 
stress and bulk modulus of the FN. This gives us the 
exact expressions 



o-fn(7) 



kfn(7) 



4/V3 kT 

(7^¥' 
4/V3 kT 



(3) 
(4) 



For the elastic constant C44, Ref. |16| finds only the 
asymptotic functional form, and therefore for our prob- 
lem we have 



C'fN44(7) 



A kT 
(7^-l)^¥' 



(5) 



with an unknown constant A. Our numerical results, pre- 
sented in Fig. ^, confirm these relations, which seem to be 
accurate over quite a large range of 7 values. In Eq. (H), 
we use the value A = 1.80 ± 0.02 obtained by fitting the 
asymptotic expression for C44 to the three data points 
corresponding to the smallest 7 values. Note that while 
in Eqs. (|^)-(@), ^hd, '«hd and Chd44 are expressed in 
units of kT/b^, in Fig. ^ they are given in units of kT/b^. 
In this representation, the stress and elastic constants of 
the FN are scaled to depend on the parameter 7 alone. 
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FIG. 2. Numerical results for the stress a (circles) , the bulk 
modulus K (squares) , and the elastic constant C44 (triangles) , 
as a function of the ratio 7 between the maximal and actual 
extensions of the net. Results are in fcT/b^ units. The solid, 
dashed and dotted curves depict the expressions on the right 
sides of Eqs. (|)-(|1), respectively [with A = 1.80 in Eq. (||)]. 
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FIG. 3. Numerical results for the stress a (circles), and the 
elastic constants Cu (squares), C12 (triangles pointing up) 
and C44 (triangles pointing down), as a function of the ra- 
tio 7 between the maximal and actual extensions of the net. 
Results are in fcT/b^ units. The lines are guides to the eye. 

Fig. m shows the dependence of the stress and elastic 
constants on 7 for weakly stretched nets. We observe 
a spectacular decay of elastic constants to almost zero 
for 7 ~ 3, and at the same time we note that the stress 
becomes independent of 7. The very fact of decrease of 
elastic constants with increasing 7 should not be surpris- 
ing, because it is intuitively clear that larger 7 represent 
a more "loose" and more "weak" solid. However, al- 



most vanishing values already at 7 '--^ 3 are not direct 
consequences of the "weakness" of the solid, but of the 
fact that a "loose" network can be approximated by a 
network of Gaussian springs. Gaussian spring is a linear 
spring of vanishing unstressed length. The energy of such 
a spring, E — ^Kr^, is simply proportional to its squared 
end-to-end distance, r^. We will show that elastic solid 
formed by such springs has exactly vanishing elastic con- 
stants, independently of the value of the spring constant 
K. Thus, the effect observed in Fig. || is an indication of 
the Gaussian nature of the system. 

It is easy to calculate the elastic properties of Gaus- 
sian networks at T = 0: The stresses of such networks 
depend on their topologies, namely on the details of the 
connectivity between the atoms and on the values of the 
springs constants between them. For 2D networks the 
stresses are not modified due to homogeneous changes in 
the size of the net, since the force exerted on the surface 
grows (diminishes) linearly with the length of the bound- 
aries. Moreover, at T = 0, the free energy, F, coincides 
with the internal energy E = Ebonds (q/j) i^a/3(r"'^)^ 
where r"^ is the length of the bond connecting atoms 
"a" and "/3" , and Ka/3 is the spring constant assigned 
to this bond. From Eq. ([|) it is obvious that the en- 
ergy expansion in the strain variables includes only lin- 
ear terms in 77, and hence, by comparing with Eq. (H), 
Cijki{T = 0) = 0. This identity, as well as the size in- 
dependence of the stresses, hold at any other tempera- 
ture since Gaussian networks have the interesting feature 
that their stress and elastic constants are temperature 
independent! For the stresses this feature is readily un- 
derstood: The stresses can be expressed as the averages 
of quantities which are linear in the coordinates of the 
atoms. When the statistical weights of the distribution 
are Gaussian, i.e., an exponent of a quadratic form of 
the coordinates, these averages coincide with the most 
probable values, namely their values at equilibrium. The 
temperature independence of the elastic constants then 
follows immediately, since the latter are just the deriva- 
tives of the stress components. 

The similarity between non-stressed tethered and 
Gaussian one-dimensional (ID) nets, i.e., linear poly- 
mers, is a consequence of central limit theorem [ pT[ . 
For topologically two-dimensional regular (non-random) 
nets, such similarity was demonstrated by Kantor et 
al. 1^: In both tethered and Gaussian two-dimensional 
nets, the mean squared distance in the embedding space, 
^xx' ~ between two distant points whose 

internal positions in the net (measured in lattice con- 
stants) are x and x', grows proportionally to In |x — x'|. 
One can define the effective spring constant, KcS^ as the 
value of if of a Gaussian network with the same connec- 
tivity and statistical properties as of the tethered net- 
work. The value of ifeff is extracted from the ratio of 
the mean squared distance, r^^, , between two points x 
and x' on the FN, and the mean squared distance f^^' 
between the same two points on a Gaussian network of 
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unit spring constants: 

^off = rJx'Axx" (6) 

f^^, can be calculated exactly, while the value of the cor- 
responding r^^, = (|r(x) — r(x')p) can be extracted from 
MC simulations of the FN with free boundaries condi- 
tions (i.e., in the absence of external pressure). We simu- 
lated a FN of 56^ = 3136 atoms and measured (using 10'' 
different configurations) r^^, for several pairs of points x 
and x' at different lattice separations. With these mea- 
surements we evaluated the effective spring constants [us- 
ing Eq. (^)], and found, as shown in Fig. ^, that for the 
FN model K^ff ~ 1.96 kT/b^. In order to support our 
conclusion about the crossover into the Gaussian regime, 
we need to show that the constant value to which the 
stress drops in Fig. |3[ is just the stress applied by a 
Gaussian net with spring constants calculated for 
non-stressed FN. For a Gaussian net with K ~ 1.96, one 
finds that a = ^/SK ^ 3.39 kT/b^ which indeed coincides 
with the value of 3.4 kT/b'^, extracted from Fig. ^. 
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FIG. 4. The effective spring constant KcB, extracted from 
MC measurements of r^^, = {|r(x) - r(x')p) [see Eq. (|)]. 
The error bars correspond to one standard deviation in the 
estimated value of r^^, . 

The persistence of the Gaussian regime to intermedi- 
ate values of 7 (7^3), is not unique for 2D nets. Such 
behavior is also found, for instance, in ID polymers. Let 
us consider, for a moment, a chain of iV ^ 1 tethers of 
maximal length 6, which is stretched by a force /, to an 
end-to-end length I = Nbg. It is a well known fact that 
this chain will be Gaussian, i.e., / and I will be propor- 
tional to each other, provided that I does not exceed the 
order of magnitude of the root mean square size of the 
chain: 

l^NbQ< y/Nb. (7) 
Yet, one must understand that in order to observe Gaus- 



sian elastic behavior, it is not essential to apply this cri- 
terion (|^) to the whole chain, but only to small segments 
of it. If there exist a certain length scale at which the po- 
tential between the atoms becomes effectively quadratic, 
i.e., can be replaced by a Gaussian spring, then the whole 
chain is like a chain of Gaussian springs, and therefore it 
is itself Gaussian. For a linear polymer chain, the effec- 
tive potential between non-neighboring atoms is calcu- 
lated by integrating out the spatial degrees of freedom of 
the atoms located between them. Such calculations are 
usually done iteratively, where on each "rescaling" step 
every second atoms is integrated out. It appears that 
even elementary potentials which are very different from 
parabola, are brought into a parabolic form within a few 
"rescaling" steps. For the specific potential used in this 
work, three steps are sufficient, which means that a seg- 
ment of iV ^ 10 tethers may be justly considered as an 
effective Gaussian spring. Similarly to a macroscopically 
large chain, we expect that the Gaussian nature of this 
segment will persist as long as it is stretched to a length 
which does not exceed its root mean square size, namely, 
as long as IO60 ^ VWb [see criterion (|^)]. This relation 
gives the lower limit, 7 = b/bo ^ y/TO ~ 3, of the Gaus- 
sian regime of a ID chain of tethers. For a 2D regular 
phantom net, the effective potential becomes approxi- 
mately parabolic also for a distance of number of bonds, 
N ^ 10 [^. Root mean square distance between two 
such points is by/In N. Thus, in order to obser ve G aus- 
sian elastic behavior, we require that IO60 < by/In 10, or, 
7 = b/bo > 10/Vln 10 ~ 4, which is consistent with the 
value 7^3, observed in Fig. ||. 

In summary, we have applied a new "fluctuation" for- 
malism to MC determination of the stress and elastic 
constants of stretched tethered networks. These systems 
provide a convenient framework for studying the entropic 
contribution to elasticity in real polymeric systems. The 
Gaussian nature of entropic elasticity, observed for non- 
stressed phantom nets, was also found when stress was 
applied. It breaks only for highly extended networks, 
close to their full-extension. This point has interest- 
ing implications to the problem of the critical elastic 
behavior of gels right above the gel-point. As already 
mentioned at the first paragraph of this paper, recently 
it was suggested by Plischke, Joos and co-workers that 
this behavior is dominated by entropy These au- 
thors studied numerically (using a different technique) 
the elastic behavior, at T ^ 0, of bond diluted (percolat- 
ing) systems at which only a fraction p of the bonds were 
present. Their results in 2D for the critical exponent / 
characterizing the growth of the shear modulus above the 
percolation threshold Pc, /i ^ {p — Pc)^ , match (within 
the range of error), the known result for the exponent t 
describing the conductivity of random resistors network 
I] (^p — The question is whether this result is uni- 
versal. For Gaussian networks the identity, f = t, can be 
proven rigorously [ p^ . One can further argue that this 
result also applies to other types of interactions, provided 
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that above a certain finite length-scale, the network is ef- 
fectively Gaussian. We have shown here that this prop- 
erty is not always insured. In a percolation problem, 
the elastic backbone is inhomogeneous and includes very 
tenuous parts where the tension applied to the network 
is distributed between very few strands. Such strands 
may deviate from Gaussian behavior when high stress is 
applied. Further complications can arise from excluded 
volume effects which have not been discussed here at all. 

This work was supported by the Israel Science Foun- 
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